{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 6.1.1 Generating surface meshes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ngsolve import *\n",
    "from netgen.csg import *\n",
    "from netgen.meshing import MeshingStep\n",
    "from ngsolve.webgui import Draw"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Closed surfaces of volume meshes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When creating a volume mesh NETGEN starts meshing the surface and afterwards the volume. By setting the perstepsend argument in the CreateMesh method we can tell NETGEN to stop after the surface is meshed resulting in a valid surface mesh."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geo = CSGeometry()\n",
    "geo.Add(Sphere(Pnt(0,0,0),1))\n",
    "mesh = Mesh(geo.GenerateMesh(maxh=0.5, perfstepsend=MeshingStep.MESHSURFACE))\n",
    "Draw(mesh)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With this arbitrary complex closed surface meshes based on CSG volume meshes can be generated."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geo = CSGeometry()\n",
    "sphere_bot = Sphere(Pnt(0,0,0),1)\n",
    "sphere_middle = Sphere(Pnt(0,1,0),0.7)\n",
    "sphere_top = Sphere(Pnt(0,1.9,0),0.4)\n",
    "geo.Add(sphere_bot+sphere_middle+sphere_top)\n",
    "\n",
    "mesh = Mesh(geo.GenerateMesh(maxh=0.5, perfstepsend=MeshingStep.MESHSURFACE))\n",
    "Draw(mesh)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Curving the resulting meshes is done as for volume meshes by"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh.Curve(2)\n",
    "Draw(mesh)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Also boundary names are specified in the same manner as for volume meshes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geo = CSGeometry()\n",
    "sphere_bot = Sphere(Pnt(0,0,0),1).bc(\"bottom\")\n",
    "sphere_middle = Sphere(Pnt(0,1,0),0.7).bc(\"middle\")\n",
    "sphere_top = Sphere(Pnt(0,1.9,0),0.4).bc(\"top\")\n",
    "geo.Add(sphere_bot+sphere_middle+sphere_top)\n",
    "mesh = Mesh(geo.GenerateMesh(maxh=0.5, perfstepsend=MeshingStep.MESHSURFACE))\n",
    "print(\"Bnd = \", mesh.GetBoundaries())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "or by calling after the mesh generation the method SetBCName"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh.ngmesh.SetBCName(0, \"other_name\")\n",
    "print(\"Bnd = \", mesh.GetBoundaries())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "where the index is 0-based. To label edges, so-called BBoundaries, if available, one can use the NameEdge method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geo = CSGeometry()\n",
    "sphere_bot = Sphere(Pnt(0,0,0),1).bc(\"bottom\")\n",
    "sphere_middle = Sphere(Pnt(0,1,0),0.7).bc(\"middle\")\n",
    "sphere_top = Sphere(Pnt(0,1.9,0),0.4).bc(\"top\")\n",
    "geo.Add(sphere_bot+sphere_middle+sphere_top)\n",
    "\n",
    "geo.NameEdge(sphere_bot,sphere_middle, \"lower\")\n",
    "geo.NameEdge(sphere_middle,sphere_top, \"upper\")\n",
    "\n",
    "mesh = Mesh(geo.GenerateMesh(maxh=0.5, perfstepsend=MeshingStep.MESHSURFACE))\n",
    "print(\"BBnd = \", mesh.GetBBoundaries())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Here, the edge resulting from intersecting two geometries is labeled. After the mesh was generated the 1-based SetCD2Name (CD = co-dimension) method can be used."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh.ngmesh.SetCD2Name(1, \"other_bbnd_name\")\n",
    "print(\"BBnd = \", mesh.GetBBoundaries())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It is possible to add and label single points, called BBBoundaries, e.g. for defining a point load, by"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geo = CSGeometry()\n",
    "geo.Add(Sphere(Pnt(0,0,0),1))\n",
    "geo.AddPoint(Pnt(0,0,1), \"pntload\")\n",
    "mesh = Mesh(geo.GenerateMesh(maxh=0.5, perfstepsend=MeshingStep.MESHSURFACE))\n",
    "Draw(mesh)\n",
    "print(\"BBBnd = \", mesh.GetBBBoundaries())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Again, renaming after mesh creation is done by (1-based index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh.ngmesh.SetCD3Name(1, \"other_bbbnd_name\")\n",
    "print(\"BBBnd = \", mesh.GetBBBoundaries())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Non-closed surfaces of volume mesh"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To generate non-closed surface meshes based on a volume geometry one has to use the AddSurface method of CSGeometry. Therefore, a (finite) volume geometry needs to be defined, which acts as base geometry. Next we use AddSurface to add a surface geometry, which is cut automatically with the base geometry. Here, no additional flag needs to be set in the GenerateMesh method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geo       = CSGeometry()\n",
    "cyl       = Cylinder(Pnt(0,0,0), Pnt(1,0,0), 1)\n",
    "bot       = Plane(Pnt(0,0,0), Vec(0,0,-1))\n",
    "right     = Plane( Pnt(3,0,0), Vec(1,0,0))\n",
    "left      = Plane(Pnt(0,0,0), Vec(-1,0,0))\n",
    "finitecyl = cyl * bot * left * right\n",
    "\n",
    "geo.AddSurface(cyl, finitecyl)\n",
    "mesh = Mesh(geo.GenerateMesh(maxh=0.5))\n",
    "Draw(mesh)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Curving, adding points, and labeling boundaries, BBoundaries, and BBBoundaries is done in the same matter as before, e.g."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geo       = CSGeometry()\n",
    "cyl       = Cylinder(Pnt(0,0,0), Pnt(1,0,0), 1)\n",
    "bot       = Plane(Pnt(0,0,0), Vec(0,0,-1))\n",
    "right     = Plane( Pnt(3,0,0), Vec(1,0,0))\n",
    "left      = Plane(Pnt(0,0,0), Vec(-1,0,0))\n",
    "finitecyl = cyl * bot * left * right\n",
    "\n",
    "geo.AddSurface(cyl, finitecyl.bc(\"surface\"))\n",
    "geo.AddPoint(Pnt(0,0,1), \"pntload\")\n",
    "\n",
    "geo.NameEdge(cyl,bot, \"sym\")\n",
    "geo.NameEdge(cyl,left, \"left\")\n",
    "geo.NameEdge(cyl,right, \"right\")\n",
    "\n",
    "mesh = Mesh(geo.GenerateMesh(maxh=0.5))\n",
    "mesh.Curve(2)\n",
    "Draw(mesh)\n",
    "print(\"Bnd   = \", mesh.GetBoundaries())\n",
    "print(\"BBnd  = \", mesh.GetBBoundaries())\n",
    "print(\"BBBnd = \", mesh.GetBBBoundaries())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Non-closed surfaces via embedding"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One can define a (possibly complicated) surfaces by a mapping\n",
    "$$ \\Phi:[0,1]\\times[0,1]\\times\\{0\\}\\to\\mathbb{R}^3$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from math import pi\n",
    "L  = 12\n",
    "W  = 1.1\n",
    "mapping = lambda x,y,z : (x*L, -W*(y-0.5)*sin(pi/2*x), W*(y-0.5)*cos(pi/2*x) )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then the MakeStructuredSurfaceMesh function is used, which takes the mapping and generates the corresponding structured surface mesh. If no mapping is given, the unit-square $[0,1]\\times[0,1]\\times\\{0\\}$ is generated."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ngsolve.meshes import MakeStructuredSurfaceMesh\n",
    "mesh = MakeStructuredSurfaceMesh(quads=True, nx=12, ny=2, mapping=mapping)\n",
    "Draw(mesh)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The number of elements in $x$ and $y$ direction (on the unit-square) can be specified, by nx and ny, as well as if a structured quadrilateral or triangular mesh should be used. The triangle orientation can be changed via the flip_triangles argument."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh = MakeStructuredSurfaceMesh(quads=False, nx=12, ny=2, mapping=mapping, flip_triangles=True)\n",
    "Draw(mesh)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note, that the boundary names are always set to \"left, right, top, bottom\" according the unit-square and can be changed via the SetCD2Name method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(mesh.GetBBoundaries())\n",
    "mesh.ngmesh.SetCD2Name(1, \"other_name\")\n",
    "print(mesh.GetBBoundaries())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To name points as BBBoundaries one can specify them in the array bbbpts and add the corresponding names in the bbbnames list. Note, that the number of elements together with the mapping have to be specified such that the points in the list are already Gridpoints. Otherwise an exception is thrown."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh = MakeStructuredSurfaceMesh(nx=12, ny=2, mapping=mapping, bbbpts=[[L,0,0]], bbbnames=[\"force\"])\n",
    "print(mesh.GetBBBoundaries())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
